Regularization parameter estimation for large-scale Tikhonov regularization using a priori information

This paper is concerned with estimating the solutions of numerically ill-posed least squares problems through Tikhonov regularization. Given a priori estimates on the covariance structure of errors in the measurement data b, and a suitable statistically-chosen σ, the Tikhonov regularized least squares functional J(σ) = ‖Ax − b‖2Wb + 1/σ 2‖D(x − x0)‖2, evaluated at its minimizer x(σ), approximately follows a χ 2 distribution with m̃ degrees of freedom. Here m̃ = m + p − n for A ∈ Rm×n, D ∈ Rp×n, matrix Wb is the inverse covariance matrix of the mean 0 normally distributed measurement errors e in b, and x0 is an estimate of the mean value of x. Using the generalized singular value decomposition of the matrix pair [W 1/2 b AD], σ can then be found such that the resulting J follows this χ2 distribution, Mead and Renaut (2008). Because the algorithm explicitly relies on the direct solution of the problem obtained using the generalized singular value decomposition it is not practical for large scale problems. Here the approach is extended for large scale problems through the use of the Newton iteration in combination with a Golub-Kahan iterative bidiagonalization of the regularized problem. The algorithm is also extended for cases in which x0 is not available, but instead a set of measurement data provides an estimate of the mean value of b. The sensitivity of the Newton algorithm to the number of steps used in the Golub-Kahan iterative bidiagonalization, and the relation between the size of the projected subproblem and σ are discussed. Experiments presented contrast the efficiency and robustness with other standard methods for finding the regularization parameter for a set of test problems and for the restoration of a relatively large real seismic signal. An application for image deblurring also validates the approach for large scale problems. We conclude that the presented approach is robust for both small and large scale discretely ill-posed least squares problems.